pacman::p_load(sf, raster, spatstat, tmap, tidyverse)Take-Home Exercise 3
Prototyping Modules for Geospatial Analytics Shiny Application
Overview
The increasing demand for eldercare facilities in Singapore highlights the importance of effective spatial planning and resource allocation. Understanding the distribution and accessibility of these facilities is crucial for both policymakers and community planners to ensure that eldercare services are adequately provided to the aging population.
This project focuses on conducting a second-order spatial point pattern analysis to investigate the spatial relationships between eldercare facilities across Singapore. By employing spatial statistical methods, we aim to identify patterns in the distribution of these facilities and assess whether they are clustered, dispersed, or randomly located within the urban landscape.
Package Installation
For this exercise, the following R packages would be used.
sf, a relatively new R package specially designed to import, manage and process vector-based geospatial data in R.
spatstat, which has a wide range of useful functions for point pattern analysis. In this hands-on exercise, it will be used to perform 1st- and 2nd-order spatial point patterns analysis and derive kernel density estimation (KDE) layer.
raster which reads, writes, manipulates, analyses and model of gridded spatial data (i.e. raster). In this hands-on exercise, it will be used to convert image output generate by spatstat into raster format.
maptools which provides a set of tools for manipulating geographic data. In this hands-on exercise, we mainly use it to convert Spatial objects into ppp format of spatstat.
tmap which provides functions for plotting cartographic quality static point patterns maps or interactive maps by using leaflet API.
Data Preparation
The datasets utilised are listed as follows:
URA’s Singapore Master Plan 2014 Subzone Boundary (No Sea) - shp format
MOH’s Location of CHAS Clinics - geojson format
MOH’s Location of eldercare services - shp format
Import Geospatial Data
Importing URA’s Singapore Master Plan 2014 Subzone Boundary (No Sea)
sg_sf <- st_read(dsn = "data/geospatial", layer="MP14_SUBZONE_NO_SEA_PL")Reading layer `MP14_SUBZONE_NO_SEA_PL' from data source
`C:\wamp64\www\crediblues\IS415-GAA\Take-Home_Ex\Take-Home_Ex03\data\geospatial'
using driver `ESRI Shapefile'
Simple feature collection with 323 features and 15 fields
Geometry type: MULTIPOLYGON
Dimension: XY
Bounding box: xmin: 2667.538 ymin: 15748.72 xmax: 56396.44 ymax: 50256.33
Projected CRS: SVY21Importing MOH’s Location of eldercare services
eldercare_sf <- st_read(dsn = "data/geospatial", layer="ELDERCARE")Reading layer `ELDERCARE' from data source
`C:\wamp64\www\crediblues\IS415-GAA\Take-Home_Ex\Take-Home_Ex03\data\geospatial'
using driver `ESRI Shapefile'
Simple feature collection with 133 features and 18 fields
Geometry type: POINT
Dimension: XY
Bounding box: xmin: 14481.92 ymin: 28218.43 xmax: 41665.14 ymax: 46804.9
Projected CRS: SVY21Importing MOH’s Location of CHAS Clinics
chasclinics_sf <- st_read("data/geospatial/CHASClinics.geojson") %>% st_transform(crs = 3414)Reading layer `CHASClinics' from data source
`C:\wamp64\www\crediblues\IS415-GAA\Take-Home_Ex\Take-Home_Ex03\data\geospatial\CHASClinics.geojson'
using driver `GeoJSON'
Simple feature collection with 1193 features and 2 fields
Geometry type: POINT
Dimension: XYZ
Bounding box: xmin: 103.5818 ymin: 1.016264 xmax: 103.9903 ymax: 1.456037
z_range: zmin: 0 zmax: 0
Geodetic CRS: WGS 84Check if all 3 geospatial datasets have the same coordinate systems.
print("sg_sf")[1] "sg_sf"st_crs(sg_sf)Coordinate Reference System:
User input: SVY21
wkt:
PROJCRS["SVY21",
BASEGEOGCRS["SVY21[WGS84]",
DATUM["World Geodetic System 1984",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]],
ID["EPSG",6326]],
PRIMEM["Greenwich",0,
ANGLEUNIT["Degree",0.0174532925199433]]],
CONVERSION["unnamed",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",1.36666666666667,
ANGLEUNIT["Degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",103.833333333333,
ANGLEUNIT["Degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",1,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",28001.642,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",38744.572,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["(E)",east,
ORDER[1],
LENGTHUNIT["metre",1,
ID["EPSG",9001]]],
AXIS["(N)",north,
ORDER[2],
LENGTHUNIT["metre",1,
ID["EPSG",9001]]]]print("eldercare_sf")[1] "eldercare_sf"st_crs(eldercare_sf)Coordinate Reference System:
User input: SVY21
wkt:
PROJCRS["SVY21",
BASEGEOGCRS["SVY21[WGS84]",
DATUM["World Geodetic System 1984",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]],
ID["EPSG",6326]],
PRIMEM["Greenwich",0,
ANGLEUNIT["Degree",0.0174532925199433]]],
CONVERSION["unnamed",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",1.36666666666667,
ANGLEUNIT["Degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",103.833333333333,
ANGLEUNIT["Degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",1,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",28001.642,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",38744.572,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["(E)",east,
ORDER[1],
LENGTHUNIT["metre",1,
ID["EPSG",9001]]],
AXIS["(N)",north,
ORDER[2],
LENGTHUNIT["metre",1,
ID["EPSG",9001]]]]print("chasclinics_sf")[1] "chasclinics_sf"st_crs(chasclinics_sf)Coordinate Reference System:
User input: EPSG:3414
wkt:
PROJCRS["SVY21 / Singapore TM",
BASEGEOGCRS["SVY21",
DATUM["SVY21",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4757]],
CONVERSION["Singapore Transverse Mercator",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",1.36666666666667,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",103.833333333333,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",1,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",28001.642,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",38744.572,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["northing (N)",north,
ORDER[1],
LENGTHUNIT["metre",1]],
AXIS["easting (E)",east,
ORDER[2],
LENGTHUNIT["metre",1]],
USAGE[
SCOPE["Cadastre, engineering survey, topographic mapping."],
AREA["Singapore - onshore and offshore."],
BBOX[1.13,103.59,1.47,104.07]],
ID["EPSG",3414]]We would need to assign the correct CRS to sg_sf and eldercare_sf and reproject chasclinics_sf from wgs to SVY21 / Singapore TM
st_set_crs() and st_transform() of sf package is used as shown below.
sg_sf <- st_set_crs(sg_sf, 3414)Warning: st_crs<- : replacing crs does not reproject data; use st_transform for
thateldercare_sf <- st_set_crs(eldercare_sf, 3414)Warning: st_crs<- : replacing crs does not reproject data; use st_transform for
thatchasclinics_sf <- st_transform(chasclinics_sf,crs = 3414)
print("sg_sf")[1] "sg_sf"st_crs(sg_sf)Coordinate Reference System:
User input: EPSG:3414
wkt:
PROJCRS["SVY21 / Singapore TM",
BASEGEOGCRS["SVY21",
DATUM["SVY21",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4757]],
CONVERSION["Singapore Transverse Mercator",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",1.36666666666667,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",103.833333333333,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",1,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",28001.642,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",38744.572,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["northing (N)",north,
ORDER[1],
LENGTHUNIT["metre",1]],
AXIS["easting (E)",east,
ORDER[2],
LENGTHUNIT["metre",1]],
USAGE[
SCOPE["Cadastre, engineering survey, topographic mapping."],
AREA["Singapore - onshore and offshore."],
BBOX[1.13,103.59,1.47,104.07]],
ID["EPSG",3414]]print("eldercare_sf")[1] "eldercare_sf"st_crs(eldercare_sf)Coordinate Reference System:
User input: EPSG:3414
wkt:
PROJCRS["SVY21 / Singapore TM",
BASEGEOGCRS["SVY21",
DATUM["SVY21",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4757]],
CONVERSION["Singapore Transverse Mercator",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",1.36666666666667,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",103.833333333333,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",1,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",28001.642,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",38744.572,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["northing (N)",north,
ORDER[1],
LENGTHUNIT["metre",1]],
AXIS["easting (E)",east,
ORDER[2],
LENGTHUNIT["metre",1]],
USAGE[
SCOPE["Cadastre, engineering survey, topographic mapping."],
AREA["Singapore - onshore and offshore."],
BBOX[1.13,103.59,1.47,104.07]],
ID["EPSG",3414]]print("chasclinics_sf")[1] "chasclinics_sf"st_crs(chasclinics_sf)Coordinate Reference System:
User input: EPSG:3414
wkt:
PROJCRS["SVY21 / Singapore TM",
BASEGEOGCRS["SVY21",
DATUM["SVY21",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4757]],
CONVERSION["Singapore Transverse Mercator",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",1.36666666666667,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",103.833333333333,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",1,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",28001.642,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",38744.572,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["northing (N)",north,
ORDER[1],
LENGTHUNIT["metre",1]],
AXIS["easting (E)",east,
ORDER[2],
LENGTHUNIT["metre",1]],
USAGE[
SCOPE["Cadastre, engineering survey, topographic mapping."],
AREA["Singapore - onshore and offshore."],
BBOX[1.13,103.59,1.47,104.07]],
ID["EPSG",3414]]Geospatial Data Wrangling
SIngapore Boundary Data
Creating owin object
sg_owin <- as.owin(sg_sf)
plot(sg_owin)
summary(sg_owin)Window: polygonal boundary
69 separate polygons (24 holes)
vertices area relative.area
polygon 1 (hole) 3 -2.21090e+00 -2.83e-09
polygon 2 285 1.61128e+06 2.06e-03
polygon 3 (hole) 3 -8.83647e-03 -1.13e-11
polygon 4 668 5.40368e+07 6.91e-02
polygon 5 44 2.26577e+03 2.90e-06
polygon 6 27 1.50315e+04 1.92e-05
polygon 7 714 1.28815e+07 1.65e-02
polygon 8 (hole) 36 -4.01660e+04 -5.14e-05
polygon 9 (hole) 317 -5.11280e+04 -6.54e-05
polygon 10 77 3.29939e+05 4.22e-04
polygon 11 30 2.80002e+04 3.58e-05
polygon 12 71 8.18750e+03 1.05e-05
polygon 13 (hole) 3 -1.68316e-04 -2.15e-13
polygon 14 (hole) 36 -7.79904e+03 -9.97e-06
polygon 15 (hole) 3 -2.18000e-06 -2.79e-15
polygon 16 (hole) 3 -6.62377e-01 -8.47e-10
polygon 17 (hole) 3 -2.09065e-03 -2.67e-12
polygon 18 91 1.49663e+04 1.91e-05
polygon 19 (hole) 26 -1.25665e+03 -1.61e-06
polygon 20 (hole) 349 -1.21433e+03 -1.55e-06
polygon 21 (hole) 20 -4.39069e+00 -5.62e-09
polygon 22 (hole) 76 -1.58325e+02 -2.02e-07
polygon 23 40 1.38607e+04 1.77e-05
polygon 24 (hole) 47 -6.00395e+03 -7.68e-06
polygon 25 (hole) 12 -8.36709e+01 -1.07e-07
polygon 26 45 2.51218e+03 3.21e-06
polygon 27 142 3.22293e+03 4.12e-06
polygon 28 148 3.10395e+03 3.97e-06
polygon 29 75 1.73526e+04 2.22e-05
polygon 30 83 5.28920e+03 6.76e-06
polygon 31 211 4.70521e+05 6.02e-04
polygon 32 106 3.04104e+03 3.89e-06
polygon 33 266 1.50631e+06 1.93e-03
polygon 34 71 5.63061e+03 7.20e-06
polygon 35 10 1.99717e+02 2.55e-07
polygon 36 478 2.06120e+06 2.64e-03
polygon 37 155 2.67502e+05 3.42e-04
polygon 38 1027 1.27782e+06 1.63e-03
polygon 39 (hole) 3 -1.16959e-03 -1.50e-12
polygon 40 65 8.42861e+04 1.08e-04
polygon 41 47 3.82087e+04 4.89e-05
polygon 42 6 4.50259e+02 5.76e-07
polygon 43 132 9.53357e+04 1.22e-04
polygon 44 (hole) 3 -3.23310e-04 -4.13e-13
polygon 45 4 2.69313e+02 3.44e-07
polygon 46 (hole) 3 -1.46474e-03 -1.87e-12
polygon 47 1045 4.44510e+06 5.68e-03
polygon 48 22 6.74651e+03 8.63e-06
polygon 49 64 3.43149e+04 4.39e-05
polygon 50 (hole) 3 -1.98390e-03 -2.54e-12
polygon 51 (hole) 4 -1.13774e-02 -1.46e-11
polygon 52 14 5.86546e+03 7.50e-06
polygon 53 95 5.96187e+04 7.62e-05
polygon 54 (hole) 4 -1.86410e-02 -2.38e-11
polygon 55 (hole) 6 -7.08892e-03 -9.07e-12
polygon 56 (hole) 3 -5.55856e-03 -7.11e-12
polygon 57 234 2.08755e+06 2.67e-03
polygon 58 10 4.90942e+02 6.28e-07
polygon 59 234 4.72886e+05 6.05e-04
polygon 60 (hole) 13 -3.91907e+02 -5.01e-07
polygon 61 15 4.03300e+04 5.16e-05
polygon 62 227 1.10308e+06 1.41e-03
polygon 63 10 6.60195e+03 8.44e-06
polygon 64 19 3.09221e+04 3.95e-05
polygon 65 145 9.61782e+05 1.23e-03
polygon 66 30 4.28933e+03 5.49e-06
polygon 67 37 1.29481e+04 1.66e-05
polygon 68 4 9.47108e+01 1.21e-07
polygon 69 14672 6.97996e+08 8.93e-01
enclosing rectangle: [2667.54, 56396.44] x [15748.72, 50256.33] units
(53730 x 34510 units)
Window area = 781945000 square units
Fraction of frame area: 0.422Eldercare facilities
Converting from sf format into spatstat’s ppp format
eldercare_ppp <- as.ppp(eldercare_sf)Warning in as.ppp.sf(eldercare_sf): only first attribute column is used for
markseldercare_pppMarked planar point pattern: 133 points
marks are numeric, of storage type 'integer'
window: rectangle = [14481.92, 41665.14] x [28218.43, 46804.9] unitsplot(eldercare_ppp)
summary(eldercare_ppp)Marked planar point pattern: 133 points
Average intensity 2.632412e-07 points per square unit
Coordinates are given to 11 decimal places
marks are numeric, of type 'integer'
Summary:
Min. 1st Qu. Median Mean 3rd Qu. Max.
1 34 67 67 100 133
Window: rectangle = [14481.92, 41665.14] x [28218.43, 46804.9] units
(27180 x 18590 units)
Window area = 505240000 square unitsChecking for duplicates
any(duplicated(eldercare_ppp))[1] FALSEmultiplicity(eldercare_ppp) [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
[38] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
[75] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
[112] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1sum(multiplicity(eldercare_ppp) > 1)[1] 0There are no duplicated point events.
Combining point events object and owin object
Let’s extract eldercare events that are located within Singapore by using the code chunk below.
eldercareSG_ppp = eldercare_ppp[sg_owin]
summary(eldercareSG_ppp)Marked planar point pattern: 133 points
Average intensity 1.700887e-07 points per square unit
Coordinates are given to 11 decimal places
marks are numeric, of type 'integer'
Summary:
Min. 1st Qu. Median Mean 3rd Qu. Max.
1 34 67 67 100 133
Window: polygonal boundary
69 separate polygons (24 holes)
vertices area relative.area
polygon 1 (hole) 3 -2.21090e+00 -2.83e-09
polygon 2 285 1.61128e+06 2.06e-03
polygon 3 (hole) 3 -8.83647e-03 -1.13e-11
polygon 4 668 5.40368e+07 6.91e-02
polygon 5 44 2.26577e+03 2.90e-06
polygon 6 27 1.50315e+04 1.92e-05
polygon 7 714 1.28815e+07 1.65e-02
polygon 8 (hole) 36 -4.01660e+04 -5.14e-05
polygon 9 (hole) 317 -5.11280e+04 -6.54e-05
polygon 10 77 3.29939e+05 4.22e-04
polygon 11 30 2.80002e+04 3.58e-05
polygon 12 71 8.18750e+03 1.05e-05
polygon 13 (hole) 3 -1.68316e-04 -2.15e-13
polygon 14 (hole) 36 -7.79904e+03 -9.97e-06
polygon 15 (hole) 3 -2.18000e-06 -2.79e-15
polygon 16 (hole) 3 -6.62377e-01 -8.47e-10
polygon 17 (hole) 3 -2.09065e-03 -2.67e-12
polygon 18 91 1.49663e+04 1.91e-05
polygon 19 (hole) 26 -1.25665e+03 -1.61e-06
polygon 20 (hole) 349 -1.21433e+03 -1.55e-06
polygon 21 (hole) 20 -4.39069e+00 -5.62e-09
polygon 22 (hole) 76 -1.58325e+02 -2.02e-07
polygon 23 40 1.38607e+04 1.77e-05
polygon 24 (hole) 47 -6.00395e+03 -7.68e-06
polygon 25 (hole) 12 -8.36709e+01 -1.07e-07
polygon 26 45 2.51218e+03 3.21e-06
polygon 27 142 3.22293e+03 4.12e-06
polygon 28 148 3.10395e+03 3.97e-06
polygon 29 75 1.73526e+04 2.22e-05
polygon 30 83 5.28920e+03 6.76e-06
polygon 31 211 4.70521e+05 6.02e-04
polygon 32 106 3.04104e+03 3.89e-06
polygon 33 266 1.50631e+06 1.93e-03
polygon 34 71 5.63061e+03 7.20e-06
polygon 35 10 1.99717e+02 2.55e-07
polygon 36 478 2.06120e+06 2.64e-03
polygon 37 155 2.67502e+05 3.42e-04
polygon 38 1027 1.27782e+06 1.63e-03
polygon 39 (hole) 3 -1.16959e-03 -1.50e-12
polygon 40 65 8.42861e+04 1.08e-04
polygon 41 47 3.82087e+04 4.89e-05
polygon 42 6 4.50259e+02 5.76e-07
polygon 43 132 9.53357e+04 1.22e-04
polygon 44 (hole) 3 -3.23310e-04 -4.13e-13
polygon 45 4 2.69313e+02 3.44e-07
polygon 46 (hole) 3 -1.46474e-03 -1.87e-12
polygon 47 1045 4.44510e+06 5.68e-03
polygon 48 22 6.74651e+03 8.63e-06
polygon 49 64 3.43149e+04 4.39e-05
polygon 50 (hole) 3 -1.98390e-03 -2.54e-12
polygon 51 (hole) 4 -1.13774e-02 -1.46e-11
polygon 52 14 5.86546e+03 7.50e-06
polygon 53 95 5.96187e+04 7.62e-05
polygon 54 (hole) 4 -1.86410e-02 -2.38e-11
polygon 55 (hole) 6 -7.08892e-03 -9.07e-12
polygon 56 (hole) 3 -5.55856e-03 -7.11e-12
polygon 57 234 2.08755e+06 2.67e-03
polygon 58 10 4.90942e+02 6.28e-07
polygon 59 234 4.72886e+05 6.05e-04
polygon 60 (hole) 13 -3.91907e+02 -5.01e-07
polygon 61 15 4.03300e+04 5.16e-05
polygon 62 227 1.10308e+06 1.41e-03
polygon 63 10 6.60195e+03 8.44e-06
polygon 64 19 3.09221e+04 3.95e-05
polygon 65 145 9.61782e+05 1.23e-03
polygon 66 30 4.28933e+03 5.49e-06
polygon 67 37 1.29481e+04 1.66e-05
polygon 68 4 9.47108e+01 1.21e-07
polygon 69 14672 6.97996e+08 8.93e-01
enclosing rectangle: [2667.54, 56396.44] x [15748.72, 50256.33] units
(53730 x 34510 units)
Window area = 781945000 square units
Fraction of frame area: 0.422plot(eldercareSG_ppp)
Chas Clinics
Converting from sf format into spatstat’s ppp format
chasclinics_ppp <- as.ppp(chasclinics_sf)Warning in as.ppp.sf(chasclinics_sf): only first attribute column is used for
markschasclinics_pppMarked planar point pattern: 1193 points
marks are of storage type 'character'
window: rectangle = [0, 45475.65] x [0, 48626.7] unitsCombining point events object and owin object
Let’s extract chas clinic events that are located within Singapore by using the code chunk below.
chasclinicsSG_ppp = chasclinics_ppp[sg_owin]
summary(chasclinicsSG_ppp)Marked planar point pattern: 1192 points
Average intensity 1.524404e-06 points per square unit
Coordinates are given to 12 decimal places
marks are of type 'character'
Summary:
Length Class Mode
1192 character character
Window: polygonal boundary
69 separate polygons (24 holes)
vertices area relative.area
polygon 1 (hole) 3 -2.21090e+00 -2.83e-09
polygon 2 285 1.61128e+06 2.06e-03
polygon 3 (hole) 3 -8.83647e-03 -1.13e-11
polygon 4 668 5.40368e+07 6.91e-02
polygon 5 44 2.26577e+03 2.90e-06
polygon 6 27 1.50315e+04 1.92e-05
polygon 7 714 1.28815e+07 1.65e-02
polygon 8 (hole) 36 -4.01660e+04 -5.14e-05
polygon 9 (hole) 317 -5.11280e+04 -6.54e-05
polygon 10 77 3.29939e+05 4.22e-04
polygon 11 30 2.80002e+04 3.58e-05
polygon 12 71 8.18750e+03 1.05e-05
polygon 13 (hole) 3 -1.68316e-04 -2.15e-13
polygon 14 (hole) 36 -7.79904e+03 -9.97e-06
polygon 15 (hole) 3 -2.18000e-06 -2.79e-15
polygon 16 (hole) 3 -6.62377e-01 -8.47e-10
polygon 17 (hole) 3 -2.09065e-03 -2.67e-12
polygon 18 91 1.49663e+04 1.91e-05
polygon 19 (hole) 26 -1.25665e+03 -1.61e-06
polygon 20 (hole) 349 -1.21433e+03 -1.55e-06
polygon 21 (hole) 20 -4.39069e+00 -5.62e-09
polygon 22 (hole) 76 -1.58325e+02 -2.02e-07
polygon 23 40 1.38607e+04 1.77e-05
polygon 24 (hole) 47 -6.00395e+03 -7.68e-06
polygon 25 (hole) 12 -8.36709e+01 -1.07e-07
polygon 26 45 2.51218e+03 3.21e-06
polygon 27 142 3.22293e+03 4.12e-06
polygon 28 148 3.10395e+03 3.97e-06
polygon 29 75 1.73526e+04 2.22e-05
polygon 30 83 5.28920e+03 6.76e-06
polygon 31 211 4.70521e+05 6.02e-04
polygon 32 106 3.04104e+03 3.89e-06
polygon 33 266 1.50631e+06 1.93e-03
polygon 34 71 5.63061e+03 7.20e-06
polygon 35 10 1.99717e+02 2.55e-07
polygon 36 478 2.06120e+06 2.64e-03
polygon 37 155 2.67502e+05 3.42e-04
polygon 38 1027 1.27782e+06 1.63e-03
polygon 39 (hole) 3 -1.16959e-03 -1.50e-12
polygon 40 65 8.42861e+04 1.08e-04
polygon 41 47 3.82087e+04 4.89e-05
polygon 42 6 4.50259e+02 5.76e-07
polygon 43 132 9.53357e+04 1.22e-04
polygon 44 (hole) 3 -3.23310e-04 -4.13e-13
polygon 45 4 2.69313e+02 3.44e-07
polygon 46 (hole) 3 -1.46474e-03 -1.87e-12
polygon 47 1045 4.44510e+06 5.68e-03
polygon 48 22 6.74651e+03 8.63e-06
polygon 49 64 3.43149e+04 4.39e-05
polygon 50 (hole) 3 -1.98390e-03 -2.54e-12
polygon 51 (hole) 4 -1.13774e-02 -1.46e-11
polygon 52 14 5.86546e+03 7.50e-06
polygon 53 95 5.96187e+04 7.62e-05
polygon 54 (hole) 4 -1.86410e-02 -2.38e-11
polygon 55 (hole) 6 -7.08892e-03 -9.07e-12
polygon 56 (hole) 3 -5.55856e-03 -7.11e-12
polygon 57 234 2.08755e+06 2.67e-03
polygon 58 10 4.90942e+02 6.28e-07
polygon 59 234 4.72886e+05 6.05e-04
polygon 60 (hole) 13 -3.91907e+02 -5.01e-07
polygon 61 15 4.03300e+04 5.16e-05
polygon 62 227 1.10308e+06 1.41e-03
polygon 63 10 6.60195e+03 8.44e-06
polygon 64 19 3.09221e+04 3.95e-05
polygon 65 145 9.61782e+05 1.23e-03
polygon 66 30 4.28933e+03 5.49e-06
polygon 67 37 1.29481e+04 1.66e-05
polygon 68 4 9.47108e+01 1.21e-07
polygon 69 14672 6.97996e+08 8.93e-01
enclosing rectangle: [2667.54, 56396.44] x [15748.72, 50256.33] units
(53730 x 34510 units)
Window area = 781945000 square units
Fraction of frame area: 0.422any(duplicated(chasclinics_ppp))[1] FALSEEDA
Visualising eldercare places
tmap_mode("view")tmap mode set to interactive viewingtmap_options(check.and.fix = TRUE)
tm_shape(sg_sf)+tm_polygons()+tm_shape(eldercare_sf)+tm_dots()Warning: The shape sg_sf is invalid (after reprojection). See sf::st_is_validVisualising chas clinics
tmap_mode("view")tmap mode set to interactive viewingtmap_options(check.and.fix = TRUE)
tm_shape(sg_sf)+tm_polygons()+tm_shape(chasclinics_sf)+tm_dots()Warning: The shape sg_sf is invalid (after reprojection). See sf::st_is_validUI Design: Shiny Storyboard
.png)
Calibration Parameters
| Parameter | Type | Filter Options |
|---|---|---|
| Spatial Unit | Single Select, dropdown | region , planning area |
| Variable | Multi-select, dropdown | All, chas clinics, eldercare services |
| Area of Interest | Multi-select, dropdown | Overall Singapore , specific region , specific planning area |
pg <- sg_sf %>%
filter(PLN_AREA_N == "BUKIT MERAH")
par(mfrow=c(2,2))
plot(pg, main = "BUKIT MERAH", max.plot = 15)
pg_owin = as.owin(pg)eldercare_pg_ppp = eldercare_ppp[pg_owin]eldercare_pg_ppp.km = rescale(eldercare_pg_ppp, 1000, "km")plot(eldercare_pg_ppp.km, main="BUKIT MERAH")
2nd Order Spatial Point Pattern Analysis
Uni-variate K-function (Eldercare Facilities Only)
Assuming user chooses to view the K-function for eldercare facilities:
spatial unit: planning area
area of interest: bukit merah
correction: Ripley
Data Wrangling before K-function plot
AOI <- sg_sf %>%
filter(PLN_AREA_N == "BUKIT MERAH")
AOI_owin = as.owin(AOI)
eldercare_AOI_ppp = eldercare_ppp[AOI_owin]# Compute the K-function with Ripley correction
K_AOI <- Kest(eldercare_AOI_ppp, correction = "Ripley")
# Plot K(d) - r with additional formatting
plot(K_AOI, . - r ~ r, ylab= "K(d) - r", xlab = "Distance d (m)",
main = "Univariate K-function for Eldercare Facilities",
legend = TRUE)
abline(h = 0, col = "gray", lty = 2) # Adds a reference line at zero
Performing Complete Spatial Randomness Test
K_AOI.csr <- envelope(eldercare_AOI_ppp, Kest, nsim = 99, rank = 1, glocal=TRUE)Generating 99 simulations of CSR ...
1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20,
21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40,
41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60,
61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80,
81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98,
99.
Done.plot(K_AOI.csr, . - r ~ r, xlab="d", ylab="K(d)-r")
UI Design: Shiny Storyboard
.png)
Calibration Parameters for Uni-variate K-function
| Parameter | Type | Filter Options |
|---|---|---|
| Spatial Unit | Single-select, dropdown | region , planning area |
| Area of Interest | Multi-select, dropdown | Overall Singapore , specific region , specific planning area |
| Variable | Multi-select, dropdown | chas clinics, eldercare services |
| Correction | Single-select, dropdown | "none", "border", "bord.modif", "isotropic", "Ripley", "translate", "translation", "rigid", "none", "good" , "best" , "all" |
| Number of Simulations | Number | Number capped at 199 due to long computation time |
Bi-variate K-function (Eldercare Facilities and Chas Clinics)
Assuming user chooses to view the K-function for eldercare facilities AND chas clinics:
spatial unit: planning area
area of interest: bukit merah
Data Wrangling before K-function plot
#AOI stands for area of interest
AOI <- sg_sf %>%
filter(PLN_AREA_N == "BUKIT MERAH")
AOI_owin = as.owin(AOI)
eldercare_sf$type <- "Eldercare"
eldercare_sf_transformed <- eldercare_sf[, c("geometry", "type")]
chasclinics_sf$type <- "CHAS Clinic"
chasclinics_sf_transformed <- chasclinics_sf[, c("geometry", "type")]
# Combine eldercare and CHAS clinics into one sf object
combined_sf <- rbind(eldercare_sf_transformed, chasclinics_sf_transformed)
# Restrict to Bukit Merah (optional: ensure both are within Bukit Merah)
combined_pg <- st_intersection(combined_sf, pg)Warning: attribute variables are assumed to be spatially constant throughout
all geometries# Convert combined dataset to ppp object
combined_ppp <- as.ppp(combined_pg)Warning in as.ppp.sf(combined_pg): only first attribute column is used for
marks# Ensure "type" is a factor and set it as marks
marks(combined_ppp) <- factor(combined_pg$type)
combined_AOI_ppp <- combined_ppp[AOI_owin]BIvariate_k_AOI <- Kcross(combined_ppp, i = "Eldercare", j = "CHAS Clinic", correction = "Ripley")
# Plot K(d) - r with additional formatting
plot(BIvariate_k_AOI, . - r ~ r, ylab= "K(d) - r", xlab = "Distance d (m)",
main = "Bi-variate K-function (Eldercare Facilities and Chas Clinics)",
legend = TRUE)
abline(h = 0, col = "gray", lty = 2) # Adds a reference line at zero
Performing Complete Spatial Randomness Test
BIvariate_k_AOI_env <- envelope(combined_AOI_ppp, Kcross, nsim = 99,
i = "Eldercare", j = "CHAS Clinic",
simulate = expression(rlabel(combined_AOI_ppp)),
savefuns = TRUE)Generating 99 simulations by evaluating expression ...
1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20,
21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40,
41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60,
61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80,
81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98,
99.
Done.plot(BIvariate_k_AOI_env, . - r ~ r,
ylab = "K(d) - r", xlab = "Distance d (m)",
main = "Bivariate K-function with Envelope (Eldercare & CHAS Clinics)",
legend = TRUE)
abline(h = 0, col = "gray", lty = 2) # Adds a reference line at zero
UI Design: Shiny Storyboard
.png)
Calibration Parameters for Bi-variate K-function
| Parameter | Type | Filter Options |
|---|---|---|
| Spatial Unit | Single-select, dropdown | region , planning area |
| Area of Interest | Multi-select, dropdown | Overall Singapore , specific region , specific planning area |
| Variable | Multi-select, dropdown | All (Eldercare and Chas Clinics) |
| Correction | Single-select, dropdown | "border", "bord.modif", "isotropic", "Ripley", "translate", "translation", "none" or "best" |
| Number of Simulations | Number | Number capped at 199 due to long computation time |